Density of states for almost diagonal random matrices 
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We study the density of states (DOS) for disordered systems whose spectral statis- 
tics can be described by a Gaussian ensemble of almost diagonal Hermitian ran- 
dom matrices. The matrices have independent random entries Hi>j with small off- 
diagonal elements: (|ifj^-| 2 ) <C (\Hn\ 2 ) ~ 1. Using the recently suggested method 
of a virial expansion in the number of interacting energy levels (Journ.Phys.A 36, 



8265), we calculate the leading correction to the Poissonian DOS in the cases of 
the Gaussian Orthogonal and Unitary Ensembles. We apply the general formula to 



the critical power-law banded random matrices and the unitary Moshe-Neuberger- 



. Shapiro model and compare DOS of these models. 
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Recently, extensive attention has been devoted to unconventional random matrix theories 
(RMTs) that interpolate between the Wigner-Dyson RMT and banded RM (BRM) with the 

£ ■ 

i ■ (almost) Poissonian level statistics and can be used as a helpful tool to explore the local- 

; 

ization transition. One of these models is the power law banded random matrix (PLBRM) 

O ■ 

theoryi^ for which the variance of the off-diagonal elements reads 
X: PLBRM: m 2 ) = \ 7Z J— rrsz ■ (1) 

It is nearly constant inside the band \i — j\ < A ~ b, and decreases as a power-law function 
(iVijl 2 ) ~ j| _2a for > A. Eq. (^) is written for periodic boundary conditions of the 

PLBRM Hamiltonian. The special case a = 1 is relevant for description of critical systems 



with multifractal eigenstate a 1 ^ 3 ! 4 ' 5 , in particular for systems at the Anderson localization- 
delocalization transition point. On the other hand, it has been conjectured^ 5 that the spectral 
statistics of critical PLBRM with large b can be mapped onto the Calogero-Sutherland 
model (CS) [2J at low temperature where instead of the spectral problem one studies the 
statistics of interacting (for the real off-diagonal elements in PLBR) or non-interacting (for 
the complex off-diagonal elements in PLBR) fermions in a parabolic confinement potential. 
The case a > 1 corresponds to the power-law localization which can be found in certain 
periodically driven quantum- mechanical systems^. If a < 1/2 the spectral statistics of 
PLBRM approaches the Wigner-Dyson universality class with (3 = 1 or 2. 

The exactly solvable model of Moshe, Neuberger and Shapiro (MNS) also incorporates 
both the Poissonian and the Wigner-Dyson level statistics^. The probability distribution of 
the Hamiltonian Tt in MNS is given by P(n) = j dU Vft(n), where 

V (H) cxexpt-TrH 2 - (^)\i([uMu,H]A) ; (2) 



,2vr6, 

the matrix U is either unitary for complex Hermitian matrices TC (the unitary MNS) or 
orthogonal for real symmetric matrices "ft (orthogonal MNS), and dU is the Haar measure. 

The spectral properties of the unitary MNS turn out to be equivalent to a system of 
noninteracting one dimensional (Id) fermions in a parabolic confinement 9 . The spectral 
statistics of the orthogonal MNS coincides^ with the statistics of Id fermions in a parabolic 
potential with the long-range attractive interaction oc (xj — Xj)~ 2 . This model of strongly 
correlated fermions is a particular case of the Calogero-Sutherland model which has been 
intensively studied as a toy model for the fractional statistics. In both cases the parameter 
b of MNS corresponds to the inverse temperature of CS: b ~ 1/Tcs- 

The connection between the two models is especially clear in the unitary case where 
the unitary matrix U = M diag{e'" Pi } can be diagonalized by a unitary transformation. 
Then the variances of V id = (m^HM) in MNS are given by: 

MNS : (IVy 2 ) = ^ . (3) 

2 l + (£) 

One can easily see that Eq.(jH} coincides with Eq.Q at a = 1 if the phases = 2im/N are 
arranged as an ordered array on a circle. In general the MNS model can be considered as 
an extension of the PLBRM model for the case of the arbitrary arrangement of phases ip n 
homogeneously distributed over the circle. 



The following formula is valid to calculate the averaged value of an observable A(H) 
which is invariant under the transformation H — > M^TiM: 

{{A){l)(J = JP(M})DM} • (4) 



Here P({<£>}) is the joint probability distribution of phases 
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and {A) ^ stands for the averaging over the Gaussian random matrix H with entries having 
zero mean value and the variance given by Eq. (JSJ • 

The two-point correlation function, which follows from Eq.(jHJ) after the integration over all 
but two phases, was calculated by Gaudin with the help of the model of free non-interacting 
fermions with a linear spectrum^: 

2 

, s = <l> i -^ = ((p i -(p j )(N/2ir). (6) 



n 2 { s ) = i - 1 



e 1 b du 



(2tt6) 2 7_iog( e 2-i>_i) + 1 

If \s\ ^> b, the correlation function is almost constant M ^ I ^ bj — > 1. There is a repulsion 
between phases at a small scale controlled by b: ^(j^l <C bj ~ (s/b) 2 . 

For b ^> 1 the spectral statistics of the critical PLBRM with a = 1 and MNS are 
asymptotically the same and at b — > oo they approach the Wigner-Dyson statistic o 2 ^ 1 ^ . 
This is because the phase repulsion in MNS is strong at large b. The phases <frij form an 
approximately equidistant lattice-like structure^. In the opposite case b <C 1 , the phase 
repulsion in MNS is weak and the phases 0^ do not form a regular structure. Disorder in 
the phase arrangements at a small distances \<f>% — <pj\ ~ 1 may become especially important. 
Therefore, there is no evident correspondence between critical PLBRM and MNS at b <C 1 
even though both ensembles have the same first correction to the Poissonian level rigidity^*^ 
and the numerics have revealed a relatively small difference in the level rigidity of PLBRM 
and MNS at b ~ 1 Q. 

The progress in BRM and PLBRM theories became possible because of mappingi*^ onto 
the nonlinear supersymmetric sigma-model-^ that allowed to obtain rigorous results by using 
various powerful methods of the field theory. However, the sigma-model always starts from 
delocalized (i.e. diffusive) modes and such mapping is only justified if the bandwidth A> 1. 
In the opposite case where all the off-diagonal matrix elements are parametrically small 



compared to the diagonal ones and the system is close to the localization, no field-theoretical 
approach is known so far. Yet such almost diagonal RMT may possess nontrivial properties 
because of the slow decay of the off-diagonal matrix elements (\Hij\ 2 ) with increasing \i— j\. 
For instance it is of fundamental interest to study the spectral statistics in systems with 
power-law localization that takes place in the power-law banded random matrix ensembles 
at a > 1. Another problem to study is the critical almost diagonal PLBRM. It is known 
that the eigenvectors of PLBRM with a = 1 remain multifractal for an arbitrary small 
value of A 5}. This means that the typical eigenfunction is extended though very sparse 
at small A. Their fractal dimensions are small as compared to the dimension d = 1 of the 
underlying chain with the long range hopping. Thus almost diagonal PLBRMs may display 
the localization- derealization transition with changing the exponent a as well as their large 
bandwidth counterpart^. 

Recently we have suggested a new method that allows to study spectral statistics of a 
disordered system described by an ensemble of the almost diagonal random matrices^. It 
is a virial expansion in the number of interacting energy levels. Unlike the field-theoretical 
approach, the virial expansion starts from the Poissonian statistics and yields a regular 
expansion in powers of the small parameter controlling the ratio of the off-diagonal elements 
to the diagonal ones (\H i9 Lj\ 2 ) / (H^) ~ b 2 <C 1. The expansion has been represented by 
the summation of diagrams which are generated with the help of the Trotter formula. A 
rigorous selection rule has been established for the diagrams, which allows to account for 
exact contributions of a given number of resonant and non-resonant interacting levels. The 
method offers a controllable way to find an answer to the question when a weak interaction 
of levels can drive the system from localization toward criticality and derealization. An 
example of the spectral form-factor has been considered for a generic dependence of the 
variance (|ifj^j| 2 ) on the difference % — j. It has been shown that a term of the order of 
foc-i j g governed by the interaction of c energy levels. The general theory has been applied 
to the Rosenzweig-Portep^ model and to the critical PLBRM. 

In the present paper, we continue studying the spectral statistics with the help of the 
virial expansion. We calculate the density of states (DOS) for the ensembles of the Gaussian 
almost diagonal random matrices. Based on the detailed presentation of the method in the 
paper Q], we will explain the corresponding diagrammatic technique for DOS. We derive 
the leading correction to the Poissonian DOS for the models of critical PLBRM and MNS. 



Let us consider a Hermitian random matrix (RM) of size NxN, N ^> 1, from a Gaussian 
ensemble. We assume that the matrix entries are random and independent. The RM is the 
Hamiltonian H of the matrix Schrodinger equation Hip n = e n ip n , where e n and ip n are 
the eigenvalues and eigenvectors, respectively. We define statistical properties of the matrix 
entries as: 

(Bu) = 0\ {Hl) = -^ (\H^) = - j\) , i ± j ; (7) 

where J-{\i — j\) > is a smooth function of its argument, and the parameter b is small: 

6< 1. 

The condition i <C 1 means that RM is almost diagonal. The parameter (3 corresponds 
to the Dyson symmetry classes: Pgoe = 1, Pgue = 2. The brackets (...) denote the 
ensemble averaging. 

We will study the spectral properties of the system concentrating on the ensemble aver- 
aged density of states: 

p(E) = (£5(E-e n )). (8) 

n 

For almost diagonal RM the representation of spectral statistics in the time domain is more 
convenient therefore we explore below the Fourier transformed DOS as a function of time: 

C(t) = (Tre iflt ). (9) 

We start with a brief explanation of the method of the virial expansion that has been 
developed in detail in Rcf.[15]. As far as we investigate the properties of almost diagonal 
RMs the Hamiltonian can be naturally divided into a diagonal part H E and a matrix of 
hopping elements V: 

H = H e + V. (10) 

For a strictly diagonal matrix (V = 0) the Poissonian DOS Cp(t) is calculated straightfor- 
wardly from Eq.Q: 

C P (t) = Ne~^ . (11) 

It follows from the definition that the hopping elements Hij = Vij ~ b are small com- 
pared to the diagonal ones H^i = £j ~ 1. However, a direct expansion of the exponentials 
e \(H e +v)t j n Eq.(j2J j n terms of V would involve serious difficulties because the matrices 



H e and V do not commute with each other. One possible way to overcome these prob- 
lems is to represent ^(He+v)* as a product of exponentials containing matrices H e and V 
separately. We do this using the the Trotter formula^: 

n 

e A+B = Um ( e A/n e B/n\ ^ ^Ht = Um TT LiH e t/n e iVt/n\ _ (12) 

p=l 

In order to obtain corrections to C(t) proportional to (bt) m one has to expand m different 
exponentials in the infinite product in the r.h.s. of Eq. (jl2J) . e ltv ^ n ~ 1 + i— , setting in the 
rest n — m exponentials V — > 0, and to perform the Gaussian averaging over H e and V. 
The terms with an odd power m give zero after the disorder averaging. Thus, the power 
m must be even and we can substitute m — > 2k. 

There are 2 k\{2k-n)\ wa Y s to choose 2k exponentials to be expanded in V from the r.h.s. 
of Eq. (jl2|) . Therefore, before taking the limit n — > oo, one has to account for all possible 
different arrangements of the expanded exponentials which results in a summation over 
the Trotter variables p\,p2, ■ ■ -P2k-\] Yh=v 1 Pi = n - Each discrete variable pi denotes the 
number of successive exponentials e lHet ^ n fused together: 

(. . . e ifl ^ n (l + if ) X e iH e t/n ^ iH e t/n ^ ^ + { t± J g i H e t/n _ ^ = 

Pi exponentials where v^o 

— e 'iH £ tpi_ 1 /n ^ _|_ |tvA e iH^tpi/n ^ _|_ g i # e p ;+:L t/n ^ 

We can introduce scaled variables Yj = p«/n converting the summation over pi to the 
integration over Y\ and eliminating the parameter n from further calculations. The resulting 
expression must be averaged over the diagonal elements Si and yields the integral I pit, k) 
which depends on the time t, on the power k and on the parameter j3. 

For almost diagonal RMs, the higher number of the interacting energy levels the smaller 
is the correction to the Poissonian level statistics governed by that interaction. In particular, 
to find the leading in b correction Ci(t) to Cp(t), 

C(t) ~ C P (t) + d(t) + ...; d(t, b < 1) < C P (t) , 

we retain in the obtained series all terms that correspond to an interaction of two different 
species of the diagonal elements Si and Ej via the hopping elements Vij and Vj^ for any 
indices i ^ j in the range from 1 to N. Then at fixed k we find the following contribution 




FIG. 1: Graphic illustration of the series Iflfijl for the leading correction C\(t) to the Poissonian 
DOS Cp(t). Shadowed boxes mark the energy levels with different shadowing ("colors") for e% 
and £j. In each diagram with a given k, they are connected by the k interaction lines which are 
associated with the factor t 2k (V 2k ) ~ (tb) 2k . 

to the correction C\(t): 

N , 

c[ k) = (it) 2k I,(t,k) £<(V^) >. 

After integration over Ik — 2 Trotter variables {Yi, Y 2 , . . . Y 2 k-2}i 1/3 can be simplified to 
one-dimensional integral: 

l,(t,k) = (1/4-^r^dF. (13) 

In accordance with the definition (J7J), the Gaussian average of ({Vi,jVj^\ ) can be trans- 
formed to the following form: 

((v i jV j ,tf)= ]C p (k)xb 2k F k (\i-j\), (14) 

where K,p{k) is the combinatorial factor. Due to the Wick theorem, it is equal to factorials: 

M*)-f (tt - 1) " ,/ '- li (») 

I*!, (3 = 2. 

(k) (k) 

We have to sum the contributions C\ over the parameter k at the end, C\ = Z^Cf ■ 
This summation yields the answer for Ci(t) as a series in powers of the product (bt): 

oo i N 

c 1 (t) = iv £(-i) fc ( bt ) 2k Mt, k) K p {k) n N (k) ; n N (k) = -J2 r\V (ie) 

k=l v£j 

see a graphic presentation in Fig^ We will show below that C\ is not larger than 0(b l ). 
We emphasize that neither the combinatorial factor Kp{k) nor the integral over the Trotter 



variables Tp(t,k) depend on the correlation function (F{\i — Thus, they are universal. 
The factor TZiy(k) is, on the contrary, model dependent. It arises because of summation 
of the product of the correlation functions — j\) over the indices i and j and is not 
universal. If we associate the Hamiltonian H with a one-dimensional chain having a long 
range hopping between different sites the summation over i and j turns out to be the 
summation in the real space along the sites of the chain. As the function (F(\i — j\) depends 
only on the difference of \i — j\ the leading part of the real space sum is: 

oo 

R N (k)~2j2F k (m) + 0(l/N). (17) 

m=l 

In what follows we will neglect the correction of the order 0(1/N) to -Rjv(fc). 

To simplify further analysis of DOS we insert Eqs. ()13ll5ll7|) into the series (fTp^l and 
change the summation order. At first, we sum over the power k getting the answer for 
Ci(t) as a one- dimensional series over m, i.e., as the sum in the real space: 

2 00 

d(t) = -2Ne~b Z(bt,m)lp(t,Z(bt,m)) , Z(bt,m) = (bt) 2 F(m) ; (18) 



1/2 



m=l 



(try , ,1 v2 



i (a-Y*)z)-h(a-Y*)z), (3 = 1 



' 7 ° I 1, (3 = 2. 



An analytical integration over Y is easily doable at (3 = 2 and the integral ()19j) simplifies 
to the error function: 



--y^^ erf (f^) . (20) 

Let us consider short and long time asymptotics. If bt -C 1, C\ is determined by a 
diagram with the minimal number of the interaction lines (see FigGJ), i-e., we can keep 
in the power series (fT7)|) the single term with k = 1 having /Cg(l) = 1, and Tp{t, 1) = 

„ ( t 2 



eX P (-2fl) erf 



2 1*| ^ V 4,3; — l/^g 



d(6t < 1) ~ -Nb 2 \t\ exp (~) erf (-^ i^(l) , (21) 



A(3) Vv4^ 

which at t < a/4/3 is parametrically smaller than Cp\ Ci/Cp|6t<i ~ & x (&£). 

One can do the Fourier transform of Eq. (j2*T|) and show that for large energies e ^> 1 
the correction to the tail of the DOS has the same Gaussian exponential dependence as the 



distribution of diagonal matrix elements of H unless -Rjv(l) is divergent: 



p(e > 1) « N\l e"^ 2 [l + b 2 f3 R N (1)] . (22) 



Thus we conclude that there is no slowly decaying Lifshitz tails for almost diagonal PLBRM 
with a > 1/2 (including the critical PLBRM with a = 1) even though the multi-point 
correlation functions may significantly deviate from the Poisson distribution. 

If t 3> 1, the integral Tp can be calculated approximately by substituting the Dirac 
5-function for the exponential in the integrand of Eq. (jl3j) : 



e^^-^W =► I p (t,k)*e-V ^r 1 2M _ lfc Y ,( fc _ 1) , +0(l/N 2 )j , (23) 

and we obtain the simplified version of the series f)18|) : 

r— _£ ~Z(6t,m) _£C^™) f / (^(6t,m))-Ji(z(6t,m)), = 1; 
Ci(t > 1) ~ -AT V/37re 4 " V — ^ — -e * x< v y v y 

m=l 1*1 [ 1, p = 2. 

(24) 

We note that ~ 6 x (6t) and, thus, we can schematically rewrite the asymptotic 

expression ()24j) as Ci(t > 1) ~ & x C\(bi). If the function Ci(bt) is finite at any value of 
(bt) the correction Ci is again parametrically smaller than Cp (see examples below). 

Now we will apply the general formulae (|lfiJ18j) for the specific models of our interest 
with the definite correlation function JF(|z — 

The model of Power law banded random matrices. The model of PLBRM is defined by 
Eq.(fTf . We restrict ourselves to the critical and almost diagonal PLBRM with a = 1 and 
b <C 1 so that the variance (JH) simplifies to the following form : 

<l^f> * I ,2 , h " , + °( &4 ) • (25) 



The term of 0(b 4 ) is not important for the correction G\. We neglect this term below. 
We define the correlation function T in Eq. (|7|l for critical almost diagonal PLBRM: 

F(\i-j\) = ~ o ,i^7. (26) 

2(f) 2 W(^| ! - J |)' 

Next we note that the sum in the real space Rn, Eq. (jT7J) , is governed by small distances 
m N and, therefore, its leading term does not depend on the boundary conditions for the 



underlying chain with the long range hopping: 



N N i oo 



E F k H = E ¥ fN ^ 2kf \ * E 7^T k + 0(l/N) . (27) 
m=i m=i 2 (fj sin 2fc (fm) m =i(2m 2 ) 



Substituting Eq. (j2*7j) into the series (fTHj) we find the correction Ci: 









) ^ 











(28) 



Let us consider the limiting cases of small and large time. If bt <C 1, we insert 
7£jv(l) = tt 2 /6 into Eq. (|2*Tj) and arrive at: 

C x {bt « 1) * -iV^^& 2 1*1 exp (-Q erf (-|=) . (29) 



We emphasize that the sum in real space 7^(1) = J2m m ~ 2 converges at m ~ 1. Thus, 
when time t is small compared to b^ 1 the correction to the Poissonian DOS of the critical 
PLBRM is sensitive to the behaviour of the correlation function T{m) at short distances. 
This statement is in fact more general and holds true for any ensemble of the almost diagonal 
PLBRM where a > 1/2 and 7^(1) converges in the thermodynamical limit iV — > oo [19]. 

If t 3> 1, we substitute the correlation function (f2T)|) into the series ()24)1 which converges at 
m ~ t b. If the product (t b) is large the correction C\ is not sensitive to the short distances 
m ~ 1 and we can replace the sum over m by the integral and find the asymptotics at the 
long time t b ^> 1: 

' 2, (3 = 1; 



t 2 



d(tb > 1) ~ -Nbe 43 x < 



(30) 

7T, = 2. 



The Mo she- Neuberger- Shapiro model. The model of MNS is defined by Eq.©. Let us 
consider the unitary case with (3 = 2. From Eqs.tjHJ) and (J7|) we find the function T for MNS: 



H\<t>i - 0,1) = \ — - — r — • (3 1 } 

2 6 2 + (f sin (f ' 



If b <C 1, the Gaudin correlation function of the phases 0jj, Eq.©, simplifies to the 
following form: 

W^iTW' (32) 

The condition (|Vij| 2 ) <C 1 holds true in MNS if |0j — 0j| > 1. However, in contrast to 
the PLBRM case where the minimal distance |0,- — 0,-1 = 1, it is violated inside the band 



< \4>i — 4>j\ < b where the matrix H is no longer almost diagonal. This band is however 
narrow at b ^ 1 . Therefore, the contribution of this band to the DOS averaged over phases 
(pi is small in the parameter b. 

We apply the strategy of the virial expansion to calculate the average over H at the 
fixed phases (pi, see Eq.fllJ). The phase averaging is done at the last step and it reduces the 
sum in real space to an integral over the difference of two phases: 

(R N (k))^ ~ 2 / F k (s) n 2 (s) ds + 0(1/N) . (33) 
Jo 

The case where the correlation function T depends on the difference of the integer indices 
can be restored from Eq.tjHHl) by substituting a sum of the 5-functions instead of the two- 
point correlator: 7^2 (s) — * Sm=i ( K S ~ m )- In full analogy with PLBRM we can prove 
that the leading term of the function {R N {k))^ in MNS does not depend on the boundary 
conditions and transform Eqs. (j31H3*3*j) to a simpler form: 

r°° f 1 1 \ k 

(R N {k)U^2j o (j^O K 2 (s)ds. (34) 

Substituting Eq. (}32|) at b <^ 1 into formula (|3*4*jl we find: 

2 b r°° 9 2 

ftl'D^wl VTsW iS - (35) 

We insert Eq. (|35jl into series and derive an analog of Eq.(|18jl for MNS: 

C l(i , = -N be -'i f (t, j^) dS , (36) 

which reduces to the following expression: 

d(t) = -TiNbe-T (l - e -T^ . (37) 

The answer (|3*7j) coincides with the leading in b term of the Fourier transformed DOS of 
MNS obtained from the model of noninteracting Id fermions in a parabolic confinement. 

Let us discuss the results obtained by means of the virial expansion for PLBRM in GUE 
(Rq. pgjl with (3 = 2) and MNS (Eq.flHEJ)). First we note that at large time scale, tb > 1, 
the summation over the real space converges at large distances (m ~ tb in PLBRM and 
s = Sb ~ tb in MNS). In this case, we can simplify the integrand in the right hand side of 
Eq. (JMj) by ignoring b in the denominator and putting 7Z 2 — > 1 • We immediately see that 



the leading term resulting from the sum over m in Eq. (j24j) is equal to the one originating 
from the integral over S in Eq. (j36|) . Thus, we conclude that the first corrections to the 
Poissonian DOS coincide for PLBRM and MNS in the long time limit^: 

tb»l =► C tf = *)\PLBRM 
Cl\MNS 

The situation is quite different in the opposite limit of the short time scale tb <C 1 where 
the asymptotics for C\ is governed by the single diagram with k = 1. This diagram is 
highly sensitive to the behaviour of the function JF at the short distances and, therefore, 
yields absolutely different answers for PLBRM and MNS after the summation in the real 
space. At a fixed value of m or s the leading diagram is of the order of ~ b 2 . In the case 
of PLBRM the sum over m reduces to a numerical prefactor in Eq.()29|) leaving the power 
of b 2 unchanged. The integration over s in the case of MNS is strongly affected by the 
region < s < b where 6 2 jF(s) TZ2(s)\ s ^b ~ 1 and the off-diagonal elements of H are of the 
order of diagonal ones. This region makes the main contribution to the correction to DOS 
which is small only because of the small volume of this region As oc b. Thus we find that 
the integration over small s < b in MNS leads to a reduction of the power of b in the short 
time limit compared to the PLBRM case: 

C\\mns 

This is a clear manifestation of the nonequivalence of PLBRM and MNS models for small 
b, or the nonequivalence of PLBRM at small b and the CS model at high temperature. 
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